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A GENERAL  TWO-DIMENSIONAL  TOKAMAK  MODEL 


INTRODUCTION 

This  report  describes  a fully  two-dimensional  computer 
simulation  model  of  tokamak  discharges.  This  Lagrangian  computer 
model  is  used  to  Investigate  the  time  evolution  of  a finite-conducti- 
vity plasma  in  a finite-conductivity  vessel,  based  on  the  quasi-static 
evolution  of  pressure  balance  In  the  plasma.  The  code  is  unique  in 
that  it  is  designed  so  as  to  follow  both  the  resistive  diffusion  and 
the  gross  dynamics  of  the  plasma. 

The  Lagrangian  coordinate  system  uses  a general-connectivity 
triangular  grid  which  permits  an  investigation  of  problems  of  much 
greater  complexity  than  conventional,  one -dimensional  schemes.  As 
examples,  non-circular  plasmas  (elliptical  or  D-shaped),  multiple 
magnetic  axes  Including  separatrices,  and  magnetic  diverters  can  be 
described  with  this  code.  For  a typical  tokamak  configuration,  the 
plasma  is  represented  by  a series  of  concentric  rings  of  vertex 
points.  These  rings  are  made  to  coincide  with  flux  surfaces  in  the 
quasi-static  approximation.  The  triangle  sides  connect  points  along 
these  rings.  More  complicated  geometries  are  just  extensions  of  this 
basic  scheme. 

Figure  1 illustrates  the  ease  with  which  a triangular  mesh 
system  can  accomodate  tokamak  features  such  as  limiters  and  con- 
ducting walls.  The  objective  of  this  model  is  to  treat  dynamics  in 
general  configurations,  including  multipoles  and  separatrices,  as  in 
Fig.  2. 

The  present  model  is  based  on  the  assumption  that  the  quasi-static 
pressure  balance  equation  > Vp  = ^l/c;  j x B)  is  always  satisfied.  The 
equations  are  solved  in  an  iterative  manner,  which  removes  the  high- 
phase-velocity  magnetosonic  and  Alfven  modes  from  consideration.  Since 
we  are  interested  in  the  long  time-averaged  behavior  of  this  transport 
model,  these  waves  are  not  of  interest. 

Note:  Manuscript  submitted  February  22,  1978. 
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The  model  incorporates  the  diffusive  fluxes  of  material,  energy 
and  magnetic  fields. 

In  the  next  section  we  present  the  mathematical  model.  The 
details  of  the  iterative  approach  to  force  balance  are  discussed  in 
the  subsequent  sections. 

MATHEMATICAL  MODEL 

The  physical  system  of  Interest  is  the  cross-section  of  a 
toroidal  conducting  shell  and  numerous  current  filaments.  The  shell 
contains  a plasma  and  vacuum  regions.  The  fields  satisfy  Maxwell's 
equations : 


V . E = 

0 

(1) 

1 ^ 

(2) 

V X E = 

c &t 

V • B = 

0 

(3) 

V X B = 

^ J 
c — 

(M 

and  J = 0 in  the  vacuum. 

Here  we  use  Gaussian  units  and  neglect  the  displacement  current, 
since  electric  current  changes  arc  assumed  slow. 

In  the  region  containing  plasma  or  conductors,  the  material 
satisfies  the  usual  magnetohydrodynamic  equations.  We  shall  assume 
the  plasma  is  everywhere  sufficiently  dense  so  that  quasineutrality 
will  hold.  Then  we  can  write  these  equations 

+ V . PV  = 0 (5) 

+ V nkT  = - J X B + M.  V (6) 

dt  — c — “ 

E + ^VxB  = TJJ  (7) 

nk  + nkV-IV  = (2-y)  nkTV*V  + hV^T  + T7j2+p,W:VV  (8) 

Where  both  viscosity  and  thermal  conductivity  will  be  important  if 
the  plasma  is  sufficiently  hot. 
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llie  equations  (1-8)  represent  the  complete  mathematical  problem 
we  would  like  to  solve.  However,  we  must  simplify  them  considerably 
in  order  to  make  the  problem  tractable  computationally. 

We  shall  assume  the  dominant  magnetic  field  arises  from  the 
current  in  the  external  conductors,  the  next  order  from  the  toroidal 
plasma  currents.  Thus  we  write  the  magnetic  field  in  terms  of  a I 

vector  potential 

B = V X A (9)  \ 


and  assume  that  A has  only  a toroidal  component,  plus  a contribu- 
tion due  to  the  toroidal  magnetic  field.  Since  the  magnetic 
potential  has  only  one  component,  we  replace  it  by  a scalar 
potential  t,  defined  by 


B 


9 9 


+ ^ X V 
2TTr  — 


(10) 


where  the  (r,  cp,  z'!  coordinate  system  has  cp  in  the  toroidal 
direction,  r the  distance  from  the  major  axis,  z the  position 
relative  to  the  midplane  of  the  system,  and  is  of  the  form 
21^/lOr  initially. 

Combining  Eqs.  2,  and  7,  the  magnetic  field  is  found  to 
satisfy 


&t 


V X (VxB)  = -^VxTJVxB. 


(11) 


After  the  use  of  several  vector  identities  and  the  approximation 
that  the  plasma  current  contributing  to  the  magnetic  field  flows  in 
the  toroidal  direction,  this  can  be  rewritten 

+ V • ( V B)  = v2  B + B • V V ( 12 ) 

This  has  the  form  of  a conservation  equation  with  diffusion.  The 
last  term  represents  the  Ware  pinch. ^ This  equation  can  also  be 
cast  in  the  form  of  a diffusion  equation  in  t,  by  substituting 
Eq.  (lO)  into  it. 
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The  complete  mathematical  statement  of  this  problem  is  now 
Eqs.  (5^  and  (12)  as  written,  and  Eqs.  (6)  and  (8)  with  J eliminated 
in  terms  of  J|.  This  is  one  conservation  equation  for  p and  three 
convection  plus  diffusion  equations  for  V,  B and  T.  They  are  to  be  i 

solved  in  a coordinate  system  determined  by  B.  The  remainder  of 

this  report  is  a description  of  the  mechanics  of  that  coordinate  I 

calculation.  J 

GEOMETRODYNAMICS 

The  magnetic  field  which  is  the  solution  of  Eq.  (12)  forms  a | 

set  of  magnetic  surfaces,  whose  locations  may  most  easily  be  deter- 
mined through  examination  of  the  values  of  >1^  in  Eq.  (lO).  By  con- 
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struction,  the  components  of  B in  the  (r,  z)  plane  lie  along  constant 
surfaces. 

For  the  problems  of  interest  to  us,  these  are  multiply- 
connected,  as  in  Fig.  2.  The  choice  of  coordinate  system  becomes 
very  important  in  such  a geometry.  Heat  and  material  flow  is  very 
rapid  along  B;  much  slower  across  it.  Thus  one  coordinate  should  be 

I 

along  the  magnetic  surfaces.  The  flux  ^ is  the  natural  coordinate 
across  them.  The  breaking  of  the  volume  into  basic  computational 
cells  is  difficult  because  of  the  presence  of  irregular 
points  and  surfaces,  such  as  the  magnetic  axes  and  separatrix. 

We  choose  to  break  the  volume  into  triangular  elements  so  that 
our  mesh  does  not  force  a spurious  structure  on  the  problem. 

We  use  methods  which  have  been  developed  previously 
for  hydrodynamics-  and  magnebohydrodynamics^.  Our  basic  cells  move 
in  the  (r,  z)  plane,  and  the  plasma  dynamics  are  calculated  by  inter- 
preting the  equations  of  motion  in  the  Lagrangian  frame  of  these 
cells.  Within  each  time  advance,  an  iteration  is  performed  to  locate 
the  cells,  then  physical  quantities  are  transported  among  the  cells. 

The  Lagrangian  pressure-balance  model  is  applied  only  to  the 
ideal  dissipation-free  pressure  and  the  J x B forces  (including  both 
the  toroidal  and  poloidal  magnetic  fields).  The  mass  in  a cell,  the 
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adiabatic  constant  and  the  magnetic  flux  in  the  cell  are  all  held 
constant  as  the  vertices  are  adjusted  in  position  in  response  to 
force  imbalance.  The  degree  of  adjustment  is  directly  proportional 
to  the  strength  of  the  force  imbalance. 

The  iterative  scheme  used  to  achieve  force  balance  is  a two- 
dimensional  Newton-Raphson  method.  The  pressure  gradient  force, 
the  X B forces  and  the  corresponding  derivatives  are  calculated 
at  each  vertex  within  the  plasma.  The  vertices  are  then  moved 
in  response  to  those  forces  while  conserving  the  mass  and  magnetic 
flux  in  a cell  and  the  adiabatic,  constant. 

Figure  5 illustrates  the  basic  structure  of  the  triangular 
grid.  If  each  vertex  is  joined  with  the  bisecting  point  on  the 
opposite  side,  then  these  lines  form  a closed  region  about  each 
vertex  as  shown.  This  region  forms  the  basic  cell  for  the  numerical 
grid.  The  pressure,  number  densities,  temperatures  and  toroidal 
field  flux  are  defined  at  these  vertices  and  assumed  constant 
throughout  the  cell.  The  toroidal  current  and  poloidal  field  flux 
are  defined  on  the  triangle  sides. 

The  V • B = 0 condition  on  the  poloidal  field  and  the  7.^=0 
condition  on  the  poloidal  current  are  identically  satisfied  by 
summing  the  corresponding  magnetic  and  current  fluxes  around  the 
triangle  sides.  In  the  vacuum  regions  the  magnetic  field  fluxes 
are  iterated  to  J = 0 while  V . B = 0. 

MAGNE  TOHYDRODYNAMI CS 

Once  the  vertex  positions  have  been  determined  to  the  required 
degree  of  convergence,  the  poloidal  and  toroidal  magnetic  diffusion 
is  determined  from  Faraday's  law  and  the  magnetic  fluxes  are  inte- 
grated forward  in  time.  The  particle  and  energy  transport  are  then 
calculated,  the  time  step  updated,  and  the  force  balance 
iteration  begins  anew. 

To  determine  the  grad  p force  at  a vertex,  we  find  the 
pressure  gradient  on  a triangle  side  common  to  two  triangles. 
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assuming  the  pressure  Is  constant  throughout  each  triangle.  (See 
Figure  4.)  The  pressure  in  a triangle  is  defined  to  be  the  average 
of  the  pressure  at  each  vertex;  i.e.^ 

(p)i  = 1/3  (pi  + p2  + Pa). 

The  volume  averaged  force  is  given  by 

1=  ! VpdV, 

and  assuming  a discrete  pressure  difference  across  some  thickness 
6 T , we  have 

F = 1/3  (P4-Pi)  2^  <r)  H 6t, 

where  6^  is  normal  to  the  triangle  side,  and  ^ is  the  length  of 
the  triangle  side.  The  (r,  z)  components  of  the  force  at  the 
vertices  at  the  ends  of  the  side  ^ are  then 

= 1/6  (P4  - Pi)  2tt  (r>  (Z2  - Z3), 

F = 1/6  (P4  - Pi)  2TT  (r)  (r2  - rg). 
z 

The  pressure  gradient  force  is  determined  and  summed  at  every  plasma 
vertex  in  a similar  manner. 

Note  that  since  the  pressure  gradient  is  calculated  as  a 
discrete  jump  across  a triangle  side  (assuming  the  pressure  is 
constant  in  each  triangle  adjoining  that  side)  if  a pressure  varia- 
tion appears  on  a flux  surface,  the  system  responds  very  slowly  to 
this  variation.  As  a result,  the  convergence  to  force  balance 
is  much  slower  along  a flux  surface  than  across  a flux 
surface. 

This  difficulty  is  removed  by  employing  a conservative  pressure- 
smoothing routine  on  the  plasma  flux  surfaces.  This  is  used  only 
when  the  plasma  has  ceased  to  compress  (or  expand);  otherwise,  the 
iteration  would  be  unstable.  The  pressure  variation  on  a flux  sur- 
face, with  the  pressure  smoothing  in  effect,  is  generally  of  the 
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order  of  1“^  and  the  algorithm  converges.  If  the  pressure  is 
smoothed  on  a flux  surface  there  is  a corresponding  adiabatic  change 
in  the  temperature. 

The  primary  poloidal  flux  variable  is  which  is  the 

poloidal  flux  through  the  i*"^  triangle  side  and  is  equal  to  the 
difference  in  the  poloidal  flux  function  between  adjacent  vertices. 
Summing  the  6\|/^  around  the  triangle  sides  assures  that  the  V • B = 0 
condition  is  identically  satisfied. 

The  magnetic  field  is  given  by  Eq.  ^lO)  and  from  Figure  5^  the 
r and  z-components  of  the  poloidal  field  for  triangle  number  one 
are  then 


B = - [(rg  - rg)  6vt/i  - (rj,  - r3)6vtfg  1 

1 (n)  Ai 

= - [(Zg  - Z3)  - (z^  - Z3)  5tg] 

^ i+TT  <ri>  Ai 


where  A^  is  the  area  of  the  triangle.  These  fields  are  assumed  to 
be  constant  throughout  the  triangle. 

The  toroidal  currents  are  determined  from  the  integral  form  of 

/ 

Ampere's  law.  The  poloidal  field  is  assumed  constant  inside  a 
triangle  and  the  toroidal  current  flowing  on  a triangle  side  is  due 
to  changes  in  the  tangential  field  across  the  boundaries.  Starting 
with  Ampere's  law  and  from  Figure  6, 

I = ^ B . d£ 

= Ett  I^B^  (ri  - rg)  + B (Zj.  - Zg ) + B Us  - ri) 

^ Z2_  LQ 

+ (Zg  - z,)l 
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and  the  current  density  is  then 


i = A ^ 


A 

In 


where  A = and  n is  the  unit  vector  in  the  toroidal  direction. 
The  poloidal  j x B force  on  a triangle  side  is  then  given  by 


F=7J’jxBdV 


= (AB  Ar  + AB  Az)  (r)  (n  x (b)  ) , 
r z “ 


where  (b)  is  the  average  poloidal  field  at  the  mid-point  of  the  side 
and  the  force  on  each  vertex  is  one-half  of  the  above  quantity. 

The  toroidal  field  (Bg)  is  stored  as  a cell-centered  flux  which 
is  independent  of  the  toroidal  coordinate  and  thus  divergence  free. 
The  poloidal  current  flux  through  the  triangle  sides  is  determined 
from  Ampere's  law 


I = ^ ^ B . di. 


= <r)  AG, 

where  A(b^)  is  the  jump  in  the  average  toroidal  field  and  the  average 
field  is  given  by 


1 3 

(B,)  - i Z 

i-l 


"9i^ 


and  thus 


A 

. In 

1 = r 


^ A(Bq>  <r>  AG 
<r>  AG  6T 


The  volume  weighted  j x B force  on  a triangle  side  due  to  the 
toroidal  field  is  then 
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F = - J j X B dV 


c — — 


2 [<r>^  Bq  ] <r>^  Bq 
^ T 


I '' 

TTT 


where 


and 


\ 


+ r.  ) 
J 


l/(^), 

i 


with 


and 


[U)^  Bq  ] = <r>^  Bg  - <r>^  Bq 

T i T1  j Tj 


<r)^Bg 


Ti  ■>  Tj 


All  the  forces  have  now  been  determined  at  the  vertices^  so  the 
vertices  can  be  moved  in  response  to  those  forces.  The  vacuum  grid 
positions  are  then  iterated  to  find  locations  giving  approximately 
equal  areas  in  the  vacuum  region.  Once  the  new  positions  have  been 
computed  throughout  the  grid,  the  new  side  lengths,  weighting  factors, 
triangle  areas  and  volumes  are  determined.  Then  the  vacuum  field 
is  iterated  to  V x ^ = 0 and  the  toroidal  field  is  made  divergence 
free  while  conserving  flux.  If  the  plasma  has  not  yet  reached  force- 
balance,  the  new  densities,  pressure  and  temperatures  are  determined 
adiabatically  (pV^=  constant)  and  subsequent  passes  are  made 
through  the  iteration  loop  until  equilibrium  is  reached. 


9 


INITIAL  CONDITIONS 


The  initial  equilibrium  position  is  determined  with  no 
toroidal  magnetic  field.  This  is  because  the  toroidal  field  is  the 
dominant  field  component  and  the  convergence  from  the  Initial  condi-  , 

tions  to  the  initial  equilibrium  state  is  an  extremely  slow  process 
in  the  full  equations.  Once  the  initial  equilibrium  state  is 

reached,  the  toroidal  and  vertical  fields  are  initialized  and  the  J 

time  evolution  begins.  | 

I 

The  iteration  normally  converges  reasonably  quickly  except  for 
low  plasma  pressure  cases  (0^  « i).  Under  these  conditions  the 
j x B forces  overwhelm  the  pressure  force  and  the  iteration  to  force 
balance  becomes  unstable. 

The  correction  to  this  problem  is  to  Include  viscosity  in  the 
basic  Iteration  to  smooth  the  vertex  motion.  The  viscosity  may  be 
viewed  as  a force  which  moves  vertices  apart  if  the  length  of  a 
triangle  side  on  a flux  surface  becomes  significantly  less  than  the 
average  length  of  the  sides  comprising  that  flux  surface.  This 
viscosity  is  very  important,  and  helps  force  the  system  to  equili- 
brium relatively  rapidly.  ' 

The  approach  to  equilibrium  proceeds  in  three  quasi-independent 
steps.  First,  the  change  in  a vertex  position  approaches  zero 
exhibiting  Vp  = 1/c  j x B force  balance  locally.  Second,  the 
average  length  of  a triangle  side  making  up  a flux  surface  approaches 
a constant  indicating  the  system  has  reached  equilibrium  globally. 

Third,  the  pressure  on  a flux  surface  approaches  a constant  (at  a 
very  slow  rate). 

Since  the  approach  to  constant  pressure  on  a flux  surface  is  a 
very  slow  process,  a conservative  pressure-smoothing  algorithm  is 
employed  to  enhance  the  approach  to  equilibrium.  This  routine  is 
not  employed  until  the  average  length  (.C)  of  the  triangle  sides 
making  up  each  plasma  flux  surface  has  reached  a constant  so  as  not 
to  destroy  convergence.  The  system  is  declared  to  have  reached  a 
steady  state  only  when  (l)  the  maximum  change  in  any  plasma  vertex 
position  ~0(10“^  cm). 
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2'  the  average  length  of  a triangle  side  on  each  plasma  flux 
surface  has  reached  a constant  value  and  the  pressure 

variation  on  a plasma  surface  is  less  than  1.5  " With  the 

toroidal  field  off,  the  system  reaches  equilibrium  after  about  100 
iterations  with  3 ~ 1 and  after  about  I'^C  iterations  with 

3 « 1.  Table  I illustrates  several  parameters  associated  with 
convergence  for  the  same  initial  conditions  with  and  without 
pressure  smoothing  after  IcO  iterations. 

TRANSPORT 

Once  the  initial  steady  state  geometry  has  been  reached,  the 
toroidal  field  and  any  vertical  fields  are  applied  and  the  physical 
d>Tiamics  begin  with  the  diffusion  of  the  magnetic  fields  accounted 
for  at  every  time  step. 

The  system  is  iterated  to  force  balance,  in  the  same  manner  as 
discussed  above,  at  every  time  step.  Once  force  balance  is  reached 
the  magnetic  fields  are  determined  in  the  plasma  and  conductors  and 
integrated  forward  in  time  with  a time  step  determined  from  the 
plasma  conductivity. 

The  magnetic  diffusion  for  the  poloidal  and  toroidal  fields  ate 
treated  separately.  In  each  case  the  electric  fields  are  determined 
algebraically  from  Ohm's  law  E = J/-.  For  the  plasma,  ^ is  the 
Spitzer  conductivity  which  is  determined  at  every  vertex  as  a 
function  of  the  temperature  and  density  at  that  vertex.  The  code 
has  provisions  to  account  for  anomalous  effects.  Once  the  electric 
fields  have  been  determined,  the  magnetic  diffusion  is  calculated 
from  the  integral  form  of  Faraday's  law. 

To  determine  the  diffusion  of  the  poloidal  magnetic  field,  the 
toroidal  electric  field  is  calculated  at  all  plasma  vertices. 

First,  the  current  and  area-weighted  conductivity  are  summed  on 
each  plasma  flux  surface 


n 
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and 


n 


o.A. 
j J 


> 


where  the  sum  is  over  the  number  of  vertices  comprising  the  flux 
surface.  The  toroidal  electric  field  at  each  vertex  is  then  given  by 

I 

r. s 


From  Faraday ' s law 

d£, 

with  the  path  of  integration  around  the  basic  cell,  the  new  poloidal 
flux  through  a triangle  side  is  given  by 
new  o Id 

6t  = - grrcAt  r.^^  - r^]. 

i-t^  i-t^ 


The  poloidal  electric  field  on  a triangle  side  is  given  by 
E = j/CT 

®i+|  “ 5tt  ^’^i'^i'^'^i  ■ ’^i+l  '^i+l^^i+l^  ^ ° 

where  is  the  toroidal  magnetic  flux  at  vertex  i.  The  new  toroidal 
flux  is  then 

new  old 

~ ~ 2T^cAt  (V  X 1^4^)* 

Once  the  magnetic  field  diffusion  has  been  calculated,  the  code 
accounts  for  particle  and  energy  transport  both  parallel  and  perpen- 
dicular to  the  magnetic  field.  (These  routines  are  still  in  the 
testing  stage.)  The  flow  velocities  miy becalculated  from  the  momentum 
equation  and  Ohm's  law,  and  the  mass  and  temperature  are  transported. 
The  model  does  not  presently  include  Internal  sinks  (a  limiter  in  the 


4 
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divertor)  so  the  transport  leaves  p and  T constant  on  flux  surfaces^ 
and  diffuses  them  across  flux  surfaces.  This  completes  the 
algorithm.  The  time  is  then  updated  and  the  new  time  step  is 
calculated  as  a fraction  of  the  resistive  diffusion  timescale.  The 
system  is  then  iterated  to  a new  equilibrium  position. 

With  the  toroidal  field  taken  into  account,  the  plasma  reaches 
force  balance  after  10-20  iterations.  More  iterations  are  required 
at  the  first  time  step. 

Both  the  approach  to  the  initial  equilibrium  state  and  the 
time  evolution  have  a built-in  check-point  restart.  Thus  the  system 
variables  are  backed  up  onto  disk  at  frequent  intervals  and  a con- 
tinual check  is  made  on  the  convergence  rate  to  equilibrium.  If, 
for  any  reason,  the  system  ceases  to  converge  properly,  the  program 
is  halted;  and  after  the  appropriate  corrections  are  made  it  can  be 
restarted.  This  also  serves  as  a protection  mechanism  for  computer 
system  crashes. 


SUMMARY  AND  CONCLUSIONS 

The  work  presented  here  is  an  outgrowth  of  the  dynamic  MHD 
Linus  2 project'^.  The  code  is  used  to  calculate  the  time  evolution 
of  a finite  conductivity  plasma  in  an  arbitrary  shaped,  finite- 
conductivity  vessel.  The  model  is  unique  in  that  it  follows  both 
the  resistive  diffusion  and  dynamical  evolution  of  the  plasma  in  an 
arbitrary  geometry.  Since  the  poloidal  electrostatic  field  is  cal- 
culated, its  role  in  the  evolution  of  trapped  electrons  can  be 
determined.  The  flexibility  of  the  triangular  grid  geometry  permits 
an  investigation  of  non-circular  tokamaks  including  Doublet  and 
divertor  systems. 

With  the  addition  of  the  check-point  restart  and  the  self- 
consistent  tests  governing  the  approach  to  equilibrium,  the  code  is 
essentially  self-sufficient.  Once  a reasonable  choice  of  initial 
conditions  is  made  (plasma  size,  geometry,  density  and  temperature; 
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toroidal  current  and  field  strength),  no  further  input  is 
required  from  the  user. 

The  introduction  of  viscosity  and  pressure  smoothing  on  a flux 

surface  admits  to  an  investigation  of  low  beta  systems  with  as  equal 

a facility  as  high  beta  (3  ~ 1)  systems  for  which  the  code  was 

originally  designed.  The  code  is  now  undergoing  thorough  testing 

by  investigating  the  outward  movement  (toward  the  limiter)  of  the 

plasma  flux  surfaces  as  3 is  increased  from  a very  small  value  to 

P 

a value  of  order  one.  These  results  are  expected  in  the  near  future. 
Future  scheduled  research  includes  general  thermal  and  particle 
transport  in  Doublet  and  divertor  systems. 
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Table  1. 


Initial  Conditions:  T = 0.1  keV,  n = I.5  x 10^^,  I = 7OKA, 


Flux  Surface  Ap(5^) 


1 

2 h.k 

5 8.1 

^ 55.5 

5 15.9 


Without  Smoothing 

(i)(cm)  maxjAr,  Azl(cm) 


2.h 

X 

10"^, 

3.7 

X 

10"® 

.660 

5.1 

X 

10"^, 

3.3 

X 

10"® 

1.6h9 

3.1 

X 

10-^, 

6.9 

X 

10"® 

1.575 

2.8 

X 

0 

1 

U1 

2.3 

X 

10"® 

2.715 

2.6 

X 

10“®, 

1.8 

X 

10"® 

With  Smoothing 

Flux  Surface  Ap(^)  (A) (cm) 


maxlAr,  Az|(cm) 


1 


2 

1.5 

.667 

3 

<1. 

1.668 

4 

1.1 

1.590 

5 

<1. 

2.731 

2.7 

X 

10"®, 

4.5 

X 

10~' 

^.3 

X 

10-®, 

2.1 

X 

10"' 

4.3 

X 

10"®, 

1.5 

X 

10"' 

5.3 

X 

10"®, 

7.7 

X 

10"' 

4.7 

X 

0 

1 

01 

4.9 

X 

10-' 
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PLASMA 

COLUMN 


Fig.  1 - Application  of  triangular  gridding  to  a typical 
tokamak  conf iguation.  Plasma  and  vacuum  regions  are  rep- 
resented as  well  as  conducting  parts  such  as  field  coils, 
conducting  wall  and  limiter. 
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Fig*  2 - Tessalation  of  the  magnetic  surfaces  in  a doublet 
configuration.  Note  the  finer  grid  near  the  two  magnetic 
axes  and  on  the  outside  in  the  regions  where  trapped-par- 
tlcle  and  ballooning  modes  are  expected. 
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Fig.  3 - Detail  of  triangular-grid  elements.  A triangle 
is  made  up  of  three  directed  line  segments,  and  a vertex 
represents  the  (shaded)  region  defined  by  the  center  of 
mass  of  each  surrounding  triangle  and  the  midpoint  of 

each  side. 
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